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Abstract 

The synchronization transition of correlated ensembles of coupled Kuramoto 
oscillators on sparse random networks is investigated. Extensive numerical 
simulations show that correlations between the native frequencies of adjacent 
oscillators on the network systematically shift the critical point as well as the 
critical exponents characterizing the transition. Negative correlations imply 
an onset of synchronization for smaller coupling, whereas positive correlations 
shift the critical coupling towards larger interaction strengths. For negatively 
correlated oscillators the transition still exhibits critical behaviour similar to 
the all-to-all coupled Kuramoto system, while positive correlations change 
the universality class of the transition depending on the correlation strength. 
Crucially, the paper demonstrates that the synchronization behaviour is not 
only determined by the coupling architecture, but is also strongly influenced 
by the oscillator placement on the coupling network. 
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1. Introduction 



The collective dynamics of synchronization has constituted a field of very 
active interest over the last couple of years. Typical areas where synchro- 
nization phenomena play an important role are manifold, ranging from fields 
as diverse as ecology, social dynamics, biological rhythms to fields as laser 
physics [H, li], Isf and also power systems 0] . These real- world systems consist 
of many elementary units that interact. Since the discovery that the cou- 
pling networks are often non-trivial jsl, loj, the research about synchronization 
phenomena on complex networks has attracted significant attention [gI, 0]- 

So far most of the work in this field has dealt with understanding the 
infiuence of the topology of the interaction network on synchronization prop- 
erties [g], 0]. Studies of this problem have mainly concentrated on three 
approaches: (i) the master stability function approach jsj to understand the 
stability of the fully synchronized state for systems of identical chaotic oscil- 
lators; (ii) several analytical methods to study the onset of synchronization 
@, [loj , and; (iii) a numerical exploration of the properties of the synchroniza- 



tion transition of network ensembles [ll|, ll2|, |13|, [IJ, |15|, |18|, ll9|, |20|. These 



studies have provided much insight, demonstating, e.g., that for symmetri- 
cal coupling small homogeneous uncliquish load-balanced networks facilitate 



the transition to complete synchronization [2l|. In asymmetrically coupled 



systems, which have mainly been investigated through weighting schemes on 



undirected networks [22|, |23|, \2% |25|, homogeneous in-signals and balanced 
loads constitute two key indicators of an enhanced synchronizability. 

In the second and third stream of this research the critical coupling that 
characterizes the onset of synchronization and its relation to structural prop- 
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erties of the coupling network have been explored, mainly using the Ku- 
ramoto model (if 



uji + a 



E 



aij sm[ 



(1) 



as a well- understood model for the study of phase synchronization (see, e.g.. 



26[ for a recent summary). In Eq. ([T]) the 0j, i = 1...N describe the phases 
of oscillators, the Ui their native frequencies, a the coupling strength and 
the matrix the top olo gy o f the coupling, i.e. the interaction network. 
Following, e.g., j§, 15, 16, 17| we set aij = 1 if i and j are connected and 
aij = otherwise. For a very good discussion of this choice of normalization 



see 



19 



For all-to-all coupling aij = 1/{N — ^ j the model ([T]) exhibits a 
second order phase transition from a desynchronized to a (partially) syn- 
chronized phase at some critical coupling strength o"c = 2/TTg{0) where 
g{-) is the distribution of the oscillator's native frequencies. Interestingly, 
the synchronization transition appears to be very similar to the mean-field 
type version of the model for some classes of complex networks such as the 
Strogatz- Watts small world model ^ij , random graphs, and even some types 



of scale- free networks [15|]. However, for some types of degree distributions of 
scale-free networks the characteristics of the synchronization transition are 
found to depend on the degree heterogeneity l3 . 

Recently, as another approach, optimization techniques were used to in- 
vestigate the relationship between correlations in oscillator placements, net- 
work architecture and synchronization 27,28,29,3^. The results from these 
studies suggest that a correlated oscillator placement can strongly affect a 
system's synchronization properties. More specifically, a grouping of oscilla- 
tors on the network in such a way that the native frequencies of linked oscil- 
lators are anti-correlated has been shown to induce a transition to complete 
synchronization for low coupling, whereas a positively correlated oscillator 
arrangement requires more coupling strength for macroscopic synchroniza- 
tion to occur. These findings have recently been corroborated by a study of 
synchronization in bipolar population networks of Kuramoto oscillators 17 . 

The results in 27, 3^, 13] suggest that the characteristics of the transition 
to synchronization may be different for positively or negatively correlated 
oscillator arrangements. However, this issue and the scaling with the system 
size has not been studied systematically so far and is the main subject of the 
present study. 



3 



In our view, these results appear particularly relevant for evolved ecolog- 
ical or biological systems, where the synchronization of the individual ele- 
ments might have been a determinant of the fitness of the system and thus 
have guided its evolution. In such a case, the evolution of the system will 
probably have been strongly affected by heterogeneity in the characteristics 
(i.e. native frequencies in the Kuramoto model) of the individual elements. 
Thus, the observed architecture of the system is not only characterized by 
the topology of the interaction network, but also strongly determined by the 
oscillator placement on the network. Similarly, as 27| highlights that gains 
in the average degree of synchronization are at least as strongly marked by 
rearrangements in the oscillator placement as by the evolution of the network 
topology itself, introducing specific oscillator correlations may be easier to 
implement than a change in the overall network arrangement in technical 
applications where synchronization between the individual elements is desir- 
able. 

Correlations between the native frequencies of adjacent oscillators can be 
measured by a standard Pearson correlation coefficient 



(2) 



where (u) = '^/Nj^i'^i is the average native frequency. It should be noted 
that strongly correlated oscillator placements are only possible on sparsely 
connected interaction networks. Generally, the more connected a network is, 
the less correlated oscillator placements are possible, i.e. for a fully coupled 



system one trivially has c^^ = 0. To our best knowledge, apart from |27|, |28 



29| . correlations between network topology and oscillator placement and their 
influence on synchronization have not found any attention in the considerable 
literature on synchronization of complex networks, cf., e.g. 0,0], for a recent 
survey. 

In this paper, we follow up on the above flndings of our previous work and 
carry out a detailed analysis of the properties of the synchronization transi- 
tion on random networks with tuneable correlated oscillator placements. As 
will be shown below, correlations in the oscillator placement can affect both 
the critical point that separates the desynchronized from the synchronized 
phase, as well as the universality class of the synchronization transition. Re- 
sults that elaborate the dependence of the properties of the transition on 
the characteristics of the oscillator arrangement on the network are detailed 
below. 
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2. Analysis of the synchronization transition 



To investigate the synchronization transition we follow a straightforward 



approach and construct random (undirected) Erdos-Renyi type graphs [31 
For technical reasons, to eliminate an irrelevant source of heterogeneity that 
otherwise would have to be averaged over, we choose the 'microcanonical' 
random graph model, where exactly L = pN distinct randomly chosen pairs 



of nodes are connected by links [32|. Starting from an initially random oscil- 
lator placement we then carry out a simple optimization procedure to gen- 
erate a correlated oscillator arrangement with a given correlation coefficient 
c* . This is achieved by randomly selecting pairs of nodes and suggesting to 
swap the associated native frequencies. Swaps are accepted if they lead to 
a correlation coefficient closer to the desired value c*, i.e. we minimize 
the difference |C(^ — c^|. The procedure is terminated if the correlation coef- 
ficient reaches the desired level of correlation within a small error tolerance, 
i.e. |ct^ — c* I <^ e, where we set e = 10"'^ for the rest of the study. Once 
the correlated oscillator ensemble is generated, we then integrate Eq.'s ([T]) 
numerically with initial conditions = 0) randomly selected from (— 7r,7r] 
and determine the standard order parameter 

r(t)e''^ = J2e"^^^'\ (3) 
j 

The amplitude order parameter r{t) is then averaged over time (after allowing 
for a sufficient time interval for relaxation) and then over several hundred 
different (correlated) oscillator Ui and network configurations aij, where the 
native frequencies Ui are selected uniformly at random from [—1, 1] [iij. More 
precisely, we measure 



i^T Kt)), (4) 



t=T„ 



where (■) indicates the average over different initial conditions and the os- 



cillator and network ensemble. As, e.g., in [ll|, Il2|, Il3|, [181, |20[ a finite size 
scaling analysis is then carried out to determine the characteristics of the 
synchronization transition. 

In the following analysis we focus on networks with average degree (k) = 
3.5, i.e. very sparsely connected networks which nevertheless already have a 
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giant component that comprises more than 95% of all nodes. Clearly, oscil- 
lators cannot be placed on loops of odd length in a perfectly anti-correlated 
manner. Thus, large densities of short loops of odd length don't allow for 
strongly anti-correlated oscillator placements on the network. In this con- 
text, the sparse connectivity ensures an asymptotically vanishing number 
of triangles and in general low densities of short loops. It is worthwhile to 
point out that in order to highlight the effect of oscillator correlations on 
the synchronization transition in a regime in which the effect is strongest, all 
experiments have been carried out on very sparsely connected networks. In 
the light of this, comparisons to previous numerical and analytical work via 



mean- field approaches [12|, Il3|, are difficult. 

Using the above optimization method correlated oscillator placements 
with Cuj in the range between —0.6 and 0.6 can be generated. Even though 
the rest of the study is focussed on only one value of network connectivity 
{k) = 3.5, we have also experimented with different link densities in this 
sparse regime. The results presented below are found to be prototypical and 
the qualitative statement appears to be generally valid for sparsely connected 
networks. 

The precise determination of the transition between the desynchronized 
and the synchronized phase requires a careful consideration of finite size 



effects. Here we follow the approach of 11, 18, i^l with the finite size scaling 
ansatz 

r(a,iV)=iV-"F((a-a,)iVi/'^), (5) 

where the exponent a = P/u is related to the exponent /3 that describes the 
scaling of the order parameter r close to the critical point in the thermody- 
namic limit, i.e. 

r ~ (a - a,f. (6) 

The exponent u in Eq. ([5]) characterizes the divergence of the correlation 
volume at the critical point and the scaling function F{x) is bounded in the 
limit X — )■ ±oo. 

From the ansatz ([5]) the critical point fXc can be determined by plots of 
rN" vs. a for varying system sizes and different values of the exponent a. 
For the correct choice of the exponent a the curves Nr"'{a) intersect at a 
unique crossing point. Having thus determined the critical coupling and the 
exponent a, the exponent l/u can be obtained as the exponent that allows 
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Figure 1: Scaling of the order parameter with the system size for various coupling strengths 
for Coj = 0.6. For a « 0.755 the decay becomes algebraic. 
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for a collapse of data points for different system sizes when plotting rN"" vs. 

Examples of this analysis are displayed in Figures [2] and [31 Thus, for 
example for Cuj = 0.6 the curves r(iV)iV" are found to intersect at a unique 
crossing point for a = 0.075 and a critical coupling of Uc = .755 is obtained 
(top left panel of Fig. This result is confirmed by plots of r{N) vs. A^, 
which shows the algebraic dependence typical for a second order transition at 
this value of a, cf. Fig. [H For larger coupling the order parameter approaches 
a value independent of the system size, whereas r{N) ~ 0{N~^^'^) for smaller 
coupling. The analysis of the data collapse close to the transition point then 
yields 1/z/ = 0.6, cf. top left panel of Fig. [31 The remainder of the transitions 
for different oscillator correlations are analysed in the other panels of Fig's [21 
and [31 

Following this FSS procedure for systematically varied values of oscilla- 
tor correlations on the networks gives the phase diagram summarized in Fig. 
[IK and the dependence of the critical exponents /3 and u on the oscillator 
correlations, cf. Fig. [3)d,c. The analysis clearly shows a monotonic de- 
pendence between oscillator correlations and critical coupling, i.e. the more 
anti-correlated the oscillator arrangement the less coupling is required for 
synchr onizat ion . 

Figure [1^ also indicates that possibilities to induce shifts to the onset 
of synchronization by oscillator anti-correlations are limited: The critical 
coupling appears to saturate for — >■ —1. Conversely, the growth of the 
critical coupling with c^^ > is well-described by an exponential dependence 
(Jcicoj) = (Tc(0)e'=-/'=o with ao = 0.96 ± .01. 

One also notes a systematic dependence of the critical exponents f3 and 
u on oscillator correlations. For c^j < our data are consistent with /3 = 0.5 
and it is evident that /3 declines monotonically with c^^ in the positively 
correlated regime. In the case of the exponent u the picture is less clear. 
Our data suggest an exponent u in the range between 2.0 < i' < 2.5 for 

< 0. Exponents for > are significantly smaller with an estimate of 

= 1.7 ± 0.1 for strongly positive correlations. 

Hence, negatively correlated and positively correlated oscillator arrange- 
ments appear to be fundamentally different regimes of the systems' organi- 
zation. For all < 0, the data indicate that the synchronization transition 
falls in the same class as the all-to-all coupled Kuramoto system, i.e. u = 2.5 
and P = 0.5 and no systematic dependence on the strength of oscillator 
anti-correlations is found. Conversely, in the regime c^j > the critical ex- 
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Figure 2: Determination of the critical points by FSS for c^j — —0.6, —0.2,0.2,0.6, panels 
(a)-(d) in the above order. The critical points obtained from the analysis are shown in 
Fig.H 



ponents vary systematically as the correlation strength is changed. Stronger 
positive correlations in the oscillator placement entail a monotonic decrease 
in P and a change in u from u = 2.2 ± 0.2 to z/ ^ 1.7 ± 0.1. Hence, for 
positively correlated oscillator placement finite size effects are found to grow 
in strength. 

2.1. Clusters of synchrony 

In this subsection we present a heuristic explanation for the observed 
dependence of the synchronization threshold on oscillator correlations. The 
explanation bases on the path towards synchronization on networks, similar 



to the findings of [19| for the evolution of clusters of synchrony on differ- 
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Figure 3: Data collapse close to the critical point for Cu: = —0.6,-0.2,0.2,0.6, panels 
(a)-(d) in the above order. The critical exponents obtained from the analysis are shown 
in Fig. m 
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Figure 4: Dependence of the (a) critical coupling, (b) critical exponent /3 and (c) critical 
exponent v on the correlations between adjacent oscillators and (d) of the exponent 7 that 
characterizes the power law scaling of relaxation times with the system size in the critical 
region. Results are from the FSS analysis, lines are merely guides for the eye. 
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ent types of networks or the observations of a dependence of the onset of 
synchronization on details of the structure of directed networks advanced in 
20 . 

We define synchronized clusters in the following way. Two oscillators 
can be considered as synchronized when the difference in their time-averaged 
frequencies cJ = limj'^^j_^oo ^/{T — T^ci) Jq <pdt is smaller than the frequency 
resolution of the numerical simulation, which is given by 1/T. Accordingly, 
we define synchronized clusters as maximum sets of synchronized pairs of 
nodes. 

Figure |5] shows simulation data for the dependence of the largest and 
the second largest synchronized cluster on the relative order parameter (a — 
ac)/(Tc- It becomes apparent that while for positively correlated oscillator 
configurations large synchronized clusters are already formed far below the 
threshold, synchronized clusters are very small below the threshold for neu- 
trally and even more so for negatively correlated configurations. Moreover, 
in strongly positively correlated ensembles, close to the transition several 
relatively large synchronized clusters not much smaller than the largest syn- 
chronized cluster exist, cf. Figures |5] and El This illustrates two different 
ways of cluster formation, similar to the differences in the formation of syn- 
chronized clusters on random graphs and scale-free networks described in 

By definition, on positively correlated oscillator configurations, oscilla- 
tors of similar native frequency tend to be neighbours. Hence, only a small 
coupling strength is required for them to synchronize. Thus, the regime be- 
low the transition to global synchrony is already marked by the presence of 
relatively large synchronized clusters. In each of these clusters oscillators 
support each other in a collective oscillation, which is mostly different from 
the global mean frequency VL = ^/N^^-cUi = 0. 

In contrast, on negatively correlated configurations the formation of syn- 
chronized clusters below the transition point is suppressed. This is due to the 
large differences in native frequencies of adjacent oscillators, which makes it 
hard to form synchronized clusters for low coupling. Further, because of the 
negative correlations of native frequencies between adjacent oscillators, each 
of the small synchronized clusters that form nevertheless, has a collective 
frequency close to the global mean = 0. 

The transition to global synchrony requires the breakup of the local clus- 
ters to form one large synchronized cluster with a collective frequency equal 
to the global mean frequency. This, however, requires more coupling strength 
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Figure 5: Dependence of the average size of the largest (LC) and second largest (SC) 
synchronized cluster on the relative order parameter (a — ac)/<Jc- The data are from 
experiments for a system of = 800 nodes and are averaged over 500 independent runs. 



for positively correlated configurations because (i) the synchronized clusters 
that have to be merged are of roughly similar size and (ii) the collective fre- 
quency of each of these clusters tends to be substantially different from the 
global mean. 

Figure [5] illustrates that this leads to two different paths to synchroniza- 
tion. In positively correlated configurations one synchronized cluster gradu- 
ally merges with smaller synchronized clusters in its vicinity. In its essentials 
this mechanism is similar to the emergence of synchronization on scale-free 



networks [19[ . On negatively correlated configurations the giant synchronized 
cluster emerges from the merger of many small synchronized clusters at the 
transition, similar to the path towards synchronization on random graphs. 

3. Relaxation dynamics 

As we pointed out before, correlations of the oscillator arrangement on 
a network also influence the timescales of the relaxation dynamics towards 
the stationary state. This problem, the relaxation dynamics on networks in 



the strongly coupled regime, has recently found attention in [3J]. Son et al. 
report that typical relaxation times scale logarithmically with the system size 
in the strongly coupled regime, but appear to follow a power law close to the 
critical point. Interestingly, this study also demonstrates that the logarithmic 
scaling is due to the heterogeneity in the initial phases and independenent of 
network characteristics. 
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Figure 6: Distribution of the sizes of synchronized clusters for {a — ac)/crc — —-05. The 
data are from experiments for a system of = 800 nodes and are averaged over 500 
independent runs. 



For our analysis here, we follow the procedure introduced in ll|, |3J| and 
renormalize the order parameter to 



r-oo - r(t) 
r = , 

where r^o = linit-^oo ''^ (^) and = r(0) and then calculate 



(7) 



''"relax 



dt'r{t' 



(8) 



As pointed out in [3J], this definition of a relaxation time has the advantage 
of being robust for obtaining a saturation time and is also useful for a system 
with many relaxation time scales. 

In Fig. [7] we plot the scaling of the relaxation times with the system 
size for two different regimes, near the critical coupling strength a ~ cxc and 
in the strongly coupled regime for a ^ Uc. Close the critical coupling, we 
observe the power law relation 



Trelax OC A^^ 



(9) 



which was already mentioned in |3J]. As can be seen from the top panel of 
Fig. [7] the exponent 7 depends on the oscillator arrangement. For strongly 
negatively correlated arrangements (c^^ = —0.6) a fit yields 7 = 0.46 ± 0.01, 
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Figure 7: Scaling of relaxation times with system size for negatively (open squares) , neu- 
trally (open circles) and positively correlated(open triangles) oscillator arrangements, (a) 
Near the critical points, power law behaviour is observed with 7 w 0.46, 0.44 and 0.66. The 
behaviour for differently correlated oscillator arrangements is distinct, (b) in the strongly 
coupled regime (a = 3. ^ cTc). The scaling of relaxation times with the system size 
is logarithmic, there is not discernable difference between relaxation times in differently 
correlated oscillator arrangements. 
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for neutral correlations (c^j = 0) a very similar exponent 7 = 0.44 ± 0.01 
is obtained, whereas for strongly positive correlations {c^ = 0.6) with 7 = 
0.66 ± .02 a clearly distinct behavior is recorded. Even though relaxation 
times are generally found to be much smaller for positively correlated oscilla- 
tor arrangements in small systems, the large scaling exponent implies longer 
relaxation times in large systems. 

The exponents for oscillator arrangements with other correlations are dis- 
played in Fig. HJi. Even though error bars are relatively large the pattern 
appears to corroborate the previous observation of the existence of two dis- 
tinct regimes of oscillator arrangements, > and c^^ < 0. For negative 
and neutral correlations the data indicate that 7 ~ 0.45 independent of the 
correlation strength, whereas a monotonic growth of 7 with c^^ is manifested 
for positive correlation strength. 

Whereas correlations in the oscillator arrangement thus clearly affect the 
behaviour in the critical region, relaxation times in the strongly coupled 
regime appear to be unaffected (bottom panel of Fig. [7]). For a ^ ac in 



agreement with the results and derivation of [3J] logarithmic scaling 



Trelax ~ ^0 + COUSt. X In (10) 

of the relaxation times with the system size is found. Data points for sys- 
tems with differently correlated oscillator arrangements almost coincide, thus 
manifesting that the relaxation behaviour far above the critical point is in- 
dependent of the oscillator correlations. 

4. Discussion 

In this paper we have discussed the synchronization phase transition of 
systems of correlated non-identical Kuramoto oscillators on very sparse ran- 
dom networks. Negative correlations between adjacent oscillators are shown 
to decrease the critical coupling and allow for a transition to partial syn- 
chronization at lower coupling strengths than positively correlated oscillator 
arrangements. The detailed analysis of the critical couplings shows that 
the benefit from negatively correlated oscillator arrangements for the onset 
of synchronization saturates quickly with the strength of anti-correlations, 
whereas the increase of the critical coupling with the positive correlation 
strength is well characterized by an exponential dependence. These results 



confirm the previous findings of [27| and demonstrate that the previous ob- 



servations are not due to the relatively small system sizes, but are valid in 
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the thermodynamic hmit. As one of the key points of the present study 
this highhghts that synchronization properties of a coupled system depend 
strongly on both, the coupling topology and the arrangment of the oscillators 
on the coupling network. 

Importantly, the analysis of the phase transitions for positively and nega- 
tively correlated oscillator arrangements reveals two fundamentally different 
regimes: (i) for anti-correlated oscillator placements the transitions is still 
characterized by the critical exponents /3 = 1/2 and v = 2.5 as for the 
uncorrected and the all-to-all coupled Kuramoto system, whereas for (ii) 
positively correlated oscillator placements the critical exponents are found 
to vary systematically with the correlation strength c^. The effect of posi- 
tively correlated oscillator arrangements is thus comparable to a change in 
the heterogeneity of the coupling network. 

We also investigated the relaxation behaviour towards the stationary 
states of the dynamics. Close to the critical point, where power law scal- 
ing of the relaxation times with the system size Tj-eiax oc A^'^ is found, oscil- 
lator correlations strongly influence the relaxation dynamics. For negative 
and neutrally correlated oscillator ensembles we find 7 ~ 0.45, whereas for 
strongly positively correlated ensembles the best estimate is 7 ~ 0.66 ±0.02. 
Conversely, far above the critical coupling strength, the scaling of relaxation 
times with the system size is logarithmic and no dependence on correlations 
in the oscillator arrangement is found. 

The results presented in this paper raise a number of interesting questions 
about the influence of oscillator correlations on the synchronization transition 
for other coupling arrangements. For instance, it would be of interest to 
investigate the influence of oscillator correlations on the critical behaviour 
on scale-free networks, a topic of quite some interest over the last years 



15L lid . I13L |18| . By the same token for the uncorrelated case it has been 
established that the lower critical dimension for phase synchronization on 
hypercubic lattices is d[ = 4, whereas it is df = 2 for frequency entrainment 



12[. Could oscillator correlations influence the lower critical dimensions for 
phase and frequency synchronization on hypercubic lattices? We think that 
both issues are interesting directions for future work. 
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